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We introduce a generalized non-uniform mean-field formalism to describe the dissociation of 
weak rod-like polyelectrolytes (PEs). Our approach allows for two-sublattice symmetry breaking 
which in titration curves is associated with a plateau for intermediate dissociation degrees. We first 
test our method in the case of a single weak PE by comparison with exact enumeration studies and 
show that it gives quantitatively accurate results for the dissociation degree in the full range of 
pH values, and in specific performs much better than the nearest-neighbor approximation (where 
exact solutions are possible). We then study charge regulation of the coupled system of a weak 
polyacid and a weak polybase as a function of their mutual distance, which has some relevance 
for PE-multilayer formation and for PE complexation. An intricate interplay of the degree of 
dissociation and the effective interaction between the PEs as a function of their mutual distance is 
found. 



I. INTRODUCTION 

The charge of weak polyacids and polybases is deter- 
mined by the probability of each functional group to dis- 
sociate and expose a charged residue. This probability 
depends on a chemical equilibrium which can be tuned 
by varying the pH of the solution. In contrast to dilute 
solutions of monoacids or monobases, in weak polyelec- 
trolytes (PEs) each functional group is influenced by all 
other groups along the polymer, via their mutual elec- 
trostatic interaction. As a result of the repulsion be- 
tween charged groups, even strong PEs become weak at 
low salt concentrations. Furthermore, when two or more 
polymers interact with each other, their degree of ion- 
ization is modified, compared to their isolated state, and 
depends on parameters such as the distance between the 
polymers and their relative spatial configuration. Due 
to the many-body nature of this problem, and the long 
range of the electrostatic interactions, an exact solution 
for the average charge as function of pH is generally not 
known. 

In this paper we consider stiff PEs, where there is no 
coupling between the dissociation degree of freedom and 
the polymer conformation (for treatment of such coupling 
in flexible PEs see, for example, We consider 

flrst a single PE in salt solution and discuss some of the 
approximations commonly used to characterize its charge 
regulation. It was previously shown that a uniform mean 
fleld approach cannot adequately describe charge regula- 
tion when the coupling between charges along the poly- 



mer is strong |4|. In these cases ionizable groups disso- 
ciate in a two-step process, characterized by a plateau 
in the charge vs. pH curve. This process results from a 
spatially inhomogeneous charging pattern and is not pre- 
dicted by a uniform mean fleld approach. We introduce 
a mean field theory with explicit symmetry breaking be- 
tween two sublattices. Such an approximation is shown 
to be semi-quantitatively accurate and performs better 
than previous calculations where the range of interac- 
tions is restricted j^], as we demonstrate by comparison 
with exact enumeration over all configurations for finite 
chain lengths. 

In the second part of the paper we apply our non- 
uniform mean-field scheme to the interaction of two stiff 
PEs. We restrict ourselves to the simple case of polymers 
aligned parallel to each other and calculate the average 
charging and free energy as function of their distance. 
Our model reveals some of the intricate effects that can 
occur in interacting weak PEs. In a broader context, 
these interactions are of interest in the formation of PE 
multilayers, composed of alternating layers of positively 
and negatively charged polymers H, IE 13 ■ particular 
weak polyacids and polybases can be used to form mul- 
tilayers la, ■ t^is '^^^ properties such as the layer 
thickness and density depend strongly on the dissocia- 
tion degree of the functional groups, and can be tuned 
sensitively by varying the pH of the solution 

II. SINGLE POLYELECTROLYTE 

The free energy for a weak PE, immersed in an aqueous 
ionic solution, can be written as follows: 
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F = -In ^ exp(-H) 

{si=0,l} 



(1) 
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The Hamiltonian Ti is given by: 



(2) 



where is the position of the ith monomer and Si can be 
either zero (for an uncharged monomer) or one (charged 
monomer). The sum in Eq. (^) goes over all different 
configurations of dissociated groups. Note that F and 
7i are given in units of the thermal energy ksT . The 
chemical potential ji is related to the pH of the solution 



I^L = -2.303(pH-pKa) -^ijK (Acid) 
^l = 2.303(pH - pKb) - Ibh (Base) 



(3) 



where k is the Debye screening length, Ib = /{ekBT) 
is the Bjerrum length, equal to about 7 A in water at 
room temperature, ksT is the thermal energy, is the 
dielectric constant of water and e is the unit charge. The 
last term in Eq. © is the self-energy of the two charges 
created in the dissociation process. We assume through- 
out this paper that the ionic solution can be described 
using the linearized Debye-Hiickel theory, so that electro- 
static interactions between charges are pairwise additive, 
as in Eq. The exact form of wdh depends on the 

salt concentration, and also on the dielectric properties 
of the polymer backbone, as will be discussed below. In 
the most simple case of dielectric continuity between the 
polymer and solution, t;DH is equal to: 



t;DH(ri,r2) = Ib- 



,-K;|ri-r2| 



(4) 



The linear Debye-Hiickel approach neglects non-linear 
effects that are associated with counterion condensa- 
tion and which are contained in the non-linear Poisson- 
Boltzmann formalism. The main reason for resorting to 
linear theory is that only at that level can the compli- 
cated problem of spatially inhomogeneous charge distri- 
butions on the PE backbone be calculated. One justifica- 
tion is that weak polyelectrolytes as studied in this paper 
are typically not strongly charged, so that non-linear ef- 
fects are less important than for strong polyelectrolytes, 
as will be discussed in more detail in the concluding sec- 
tion. 

For the following calculations, it is convenient to use 
symmetric variables Si having the values —1,1 instead of 
0,1: 



1 + Si 



In terms of these variables the partition fimction is: 



(5) 




- ^ SiSj-DDH [a(i - i) 



(6) 



where: 



j>0 



(7) 



and a is the nearest-neighbor distance between dissocia- 
ble groups on a straight line. 



A. Non-uniform mean-field approacli with two 
sublattices 

1. Mean-field equations 

In principle, the above statistical one-dimensional 
problem can be solved using transfer-matrix techniques 
which take the long-ranged interactions into account via 
a multiple-time integration with a suitably chosen kernel. 
In order to obtain a simple, tractable solution we use 
mean-field methods, which are implemented in the fol- 
lowing way. The Gibbs variational principle can be used 
to obtain an upper bound for the free energy F = — InZ, 



F<Fo + - (Ho)c 



(8) 



In this inequality TLq is a trial Hamiltonian (to be spec- 
ified below) and Fq — — ln.Zo, where Zq is the parti- 
tion function obtained from Tig; The thermal averages 
in Eq. JHJ are evaluated using Tig- We introduce the trial 
Hamiltonian 



^0 = ^0 ^ S2i + ^1 ^ S2i+1 



(9) 



which separates the polymer into two sublattices. The 
variational parameters /ig, hi are fields which act on the 
charges in the two sublattices. By minimizing the right 
hand side of Eq. © with respect to /ig and hi, the fol- 
lowing equations are obtained. 



ho = jl + J (sg)o + K (si)o 
hi = jl + J (si)o + K (sg)o 



(10) 



where 



(so)g = -tanh(ft,o) 



-tanh(/ii), (11) 



and 



1 



j>0 



vM'^ja] ; i^=i^WDH[(2j + l)a] (12) 



j>0 



The choice of two sublattices (as opposed to more sub- 
lattices with a larger period) is related to the strong anti- 
correlation that can exist between adjacent monomers, 
and will be further motivated below. 



3 



2. Main properties of mean-field equations 

Equations H10() and always have a symmetric solu- 
tion for which ho — hi. However, the symmetric solution 
is not always the minimum of the free energy but can be, 
instead, a saddle point. In these cases two other solu- 
tions exist, both of which break the symmetry between 
the two sublattices, i.e., Hq ^ hi. One solution can be 
obtained from the other by exchanging ho and hi. The 
average charging degree of the polymer is then equal to: 



(•s"o)o + (Si>o + 2 



(13) 



In order to understand for which parameters symmetry 
breaking occurs, let us consider first the case p, — 0. In 
this case the Hamiltonian exhibits the symmetry Si — s- 
—Si in addition to the symmetry of exchanging the two 
sublattices. Even if the latter symmetry is broken, we 
have (s)q = 0, or equivalently (s)q — 1/2, i.e., exactly 
half of the monomers are dissociated. Using the fact that 
ho = —hi, Eqs. H1U|) and (|ll|l reduce in this case to one 
transcendental equation, 



ho = {K - J) tanh(/io) 



(14) 



This equation has a non-zero solution (where ho ^ hi) 
only if: 



K- J >l 



(15) 



If this condition is met, a sublattice symmetry breaking 
solution also exists within a certain range of ji values 
around zero. Outside the range jl = ±flc there is no 
symmetry breaking, i.e., (so)q = (si)o- If condition 115() 
is not met, there is no symmetry breaking solution for 
any value of (1. 

The solution with ho = hi (no sublattice symmetry 
breaking) can be found by substituting this equality in 
Eas. HlU|l and l|ll|) . leading to the transcendental equation 



ho = {K + J) tanh(/io) 



(16) 



In a uniform mean-field approximation this solution is 
found for all values of fl, whereas in our case it applies 
only for \fl\ > jlc- 

Before considering concrete examples we comment on 
the nature of the transition at jl = ±flc. This transi- 
tion is second order, i.e., (si)o ~ (so)o — > as /i — > ±fj.c- 
Note that a non-zero jl does not break the symmetry of 
exchanging the two sublattices in the Hamiltonian. This 
is why the order parameter is continuous at the tran- 
sition. However the derivative of the order parameter 
with respect to p, is discontinuous and diverges when ap- 
proaching the transition from the side where symmetry 
breaking occurs. Similarly, the derivative of the average 
dissociation degree with respect to fl has a discontinuity 
at the transition. These are artifacts of the mean-field 



approach, since the exact solution for a one dimensional 
system with short-ran ged interactions cannot exhibit a 
real phase transition |ll|. However, our non- uniform 
mean-field scheme still predicts the average charge very 
accurately, as will be demonstrated below. 



B. Uniform dielectric constant 

In the case of a uniform dielectric constant, in which 
the screened interaction is given by Eq. I^J, the summa- 
tions in Eq. H12|l can be performed explicitly, yielding 



J 



Ib 
'4a 



[ln(l-e~'"°) -f ln(l + e-'"')] 



7^ = --2. [in(l-e-'=") -ln(l + e-''")] (17) 



(18) 



The condition for sublattice symmetry breaking, Eq. 
(|15|l . translates to 



K - J ^—\n(l + e~'"') > 1 
2a 



(19) 



Increasing k decreases K — J and thus the possibility for 
sublattice symmetry breaking. Symmetry breaking can 
occur for a finite range of k and pH only if the condition 



Ib J_ 
a ln2 



2.9 



(20) 



is met. For Ib = 7.0 A this condition leads to o < 2.4 A. 
Hence vinyl-based polymers with an acid group on every 
second Carbon atom such as poly-styrene-sulfonate or 
poly-acrylic-acid with a charge-distance of a ~ 2.5 A are 
marginally close to symmetry breaking within the present 
model. However, a dielectric discontinuity due to the 
polymer backbone can increase K — J considerably, as 
will be discussed in the following sub-section. 

Examples with symmetry breaking will be shown in the 
following sub-section, while here we restrict ourselves to 
the case of dielectric continuity and no symmetry break- 
ing. In the case of no symmetry breaking all the depen- 
dence on K enters through the quantity 

if + J^-(ZB/2a)ln(l-e"'"^), (21) 

which increases with increasing Debye length 

1. Results 

In Fig. 1(a) we show the average degree of charge disso- 
ciation following from our mean-field equation for a poly- 
acid with a — 2.5 A, for four different salt concentrations, 
corresponding to = 3, 10, 30 and 100 A (solid lines). 
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FIG. 1: (a) Average degree of dissociation of a polyacid as function of pH-pKa, calculated using a mean-field approximation (solid 
lines). For comparison an exact enumeration over all configurations is also shown (symbols) with A'^ = 20 and using periodic 
boundary conditions. Results for four different values of the Debye length are shown, — SA (squares), lOA (diamonds), 
30 A (triangles) and 100 A (circles). The separation between charged groups is a = 2.5 A. (b) Average degree of dissociation as 
function of pH-pKa, calculated using an exact enumeration. The solid lines show enumeration results with periodic boundary 
conditions (as shown using symbols in part (a)). These results are compared with enumeration without periodic boundary 
conditions (dashed lines). Two different values of the Debye length are shown in the plot, = 3 (to the left) and 100 A(to 
the right). 



These results are compared with an exact calculation of 
the free energy for a finite chain with = 20 dissociable 
groups (symbols), by enumeration over all 2^ states. For 
the exact enumeration, periodic boundary conditions are 
imposed by setting the interaction between monomers i 
and j to be 

oo 

<Hihj)= J2 yMhj + nN). (22) 

n— — oo 

The comparison between mean-field and the exact enu- 
meration is very good for all four values of k^^ shown in 
the figure. Note that in all these cases there is no sym- 
metry breaking in the mean-field solution, as expected 
since the charge distance of a — 2.5 A does not satisfy 
the condition Eq. (|2()|l . Comparison of the four curves 
shows that k^^ has a large effect on the degree of charg- 
ing. As is increased each monomer interacts more 
strongly with the other monomers. This increased re- 
pulsion reduces the charging, or, as one might put it, 
the long-ranged repulsion between charged groups makes 
even strong PEs weak. 

Throughout this work we will assume that the polymer 
is long compared to the Debye length. However, it is im- 
portant to realize that for shorter polymers the average 
degree of dissociation depends on the polymer length. In 
order to demonstrate this point wc show in Fig. 1 (b) the 
average degree of dissociation as function of pH for a fi- 
nite chain of length Na, where A^ = 20 and a = 2.5 A. In 
the calculation an exact enumeration over all configura- 
tions is performed, without periodic boundary conditions 
(dashed lines). The results are compared with an enu- 



meration with periodic boundary conditions (solid lines) , 
as was done in Fig. 1(a). For — SA, the polymer 
length is larger than the screening length, Na 3> 
and the two calculations yield nearly identical results. In 
the second case shown in the plot, = 100 A, and 
Na are of the same order of magnitude and there are 
significant finite size effects. These results may be im- 
portant for the interaction of short DNA oligomers with 
substrates, as they show that the polymer length affects 
adsorption behavior also via the effective charge of PEs. 



2. Restriction to nearest-neighbor interactions 

A common approximation that was previously applied 
for the charge regulation of PEs is to consider only 
nearest- neighbor (NN) interactions 0, 0, The 
reason is that exact closed-form solutions are available 
in this case. Within the non-uniform mean-field ap- 
proach, this approximation corresponds to setting K — 
{Is/ia) exp(— Ko) and J = 0. Note that the combination 
K + J, which determines the solution without symmetry 
breaking, is smaller in the NN case than in full interac- 
tion case. On the other hand the combination K — J, 
which affects the symmetry breaking (see Eq. (I15f) ). is 
larger in the NN case. 

In Fig. 2 we compare the NN predictions (dashed lines 
and open symbols) with those obtained using the full 
long-range interaction (solid lines and filled symbols), 
for two different values of the Debye length. The sym- 
bols show exact enumeration results, obtained using pe- 
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pH - pKa 

FIG. 2: Comparison of a nearest neighbor mean-field approx- 
imation (dasfied lines) with a mean-field calculation taking 
into account all long range interactions (solid lines) for a poly- 
acid having a = 2.5 A, as function of pH-pKa. Results for two 
values of the Debye length are shown. They can be distin- 
guished in the plot according to the type of symbols (circles, 
= 10 A; squares, k,~^ — 100 A). The symbols show exact 
enumeration results with A'^ = 20 and periodic boundary con- 
ditions (empty symbols - only nearest neighbor interactions; 
full symbols - full long range interaction). 

riodic boundary conditions and are thus representative 
of an infinitely long system, while the lines show mean- 
field results. For both values of the Debye length there 
are significant deviations between the NN result and the 
full interaction. As can be expected, these deviations 
are larger for the larger screening length, — 100 A, 
since in this case the interaction between further-nearest 
neighbors contributes significantly to the total interac- 
tion. Note that for = 100 A and NN interactions 
(open square symbols and dashed line) there is also a 
small effect of symmetry breaking in the mean-field solu- 
tion, resulting in a non-monotonous slope near (s) = 0.5. 
A similar effect is seen in the exact enumeration. As a 
main result of this section, we see that the restriction to 
nearest-neighbor interactions is in general not a good ap- 
proximation, while the mean-field approach reproduces 
the exact enumeration results very accurately. 

C. Non-uniform dielectric constant 

So far we neglected effects of the dielectric discontinu- 
ity between the polymer and its surroundings. As these 



effects tend to increase the electrostatic interactions, they 
are expected to be important for the dissociation process. 
They were estimated in Ref. [Jj using a simple model, 
shown in Fig. 3. The PE is modeled as a cylinder of ra- 
dius d and dielectric constant Sd < £w, where — 80 
is the dielectric constant of water. Charged groups are 

d 
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FIG. 3: Schematic representation of a simple model taking 
into account the difference between the dielectric properties 
of water and a polymer's backbone. A PE is modeled as a 
cylinder of radius d with a dielectric constant Sd, while the di- 
electric constant outside the cylinder is equal to Charged 
groups that can dissociate from the polymer are located at 
regular intervals a from each other, at a distance b from the 
cylinder axis. 



assumed to be equally spaced along the cylinder axis, 
with separation a. Here we generalize this model to 
some extent by placing these charged groups at a dis- 
tance < b < d from the axis, as shown in Fig. 3. The 
electrostatic potential exerted by one such charge on an- 
other one was calculated in Ref. .4] and is given by 

n— — oo 

where z is the distance between the charges. The first 
term is equal to the electrostatic interaction within a 
medium of dielectric constant Ed, with no screening by 
salt. In the second term, Wp is equal to 

Wp{z, n) = A—Ib / dfc cos(fcz) [In{kh)f R{k, n) (24) 

Sd Jo 

where 



R{k,n) = 



ksd [Kn-i{kd) + Kn+i{kd)] K^ipd) - pe^ [Kn-i{pd) + Kn+i{pd)] Knjkd) 
ksd [In-i{kd) + In+i{kd)] Kn{pd) +pew [Kn~i{pd) + Kn+i{pd)] I„{kd) 



(25) 
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1/2 

p = {k^ + K^) and Km In are the n-th modified Bessel 
functions of the first and second kind, respectively. 

In Ref. jj] only the case 6 = was considered. For 
two charges located exactly at the polymer axis, 6 = 0, 
the electrostatic interaction approaches {Sw / £d)lB / z at 
short separations, while for larger separations it crosses 
over to the interaction in the aqueous ionic solution, 
ZBexp(— The former interaction is typically much 
larger than the latter, leading to a two-step charging 
curve and failure of a uniform mean-field approach. On 
the other hand, only close-by monomers interact strongly 
with each other, motivating the use of a nearest-neighbor 
(NN) model where only interactions between neigh- 
boring monomers are taken into account. The free en- 
ergy and avera ge d egree of dissociation can then be cal- 
culated exactly JT^ IT^ . e.g., using the transfer matrix 
method [Ulla. 

The NN approximation indeed predicts the two-step 
behavior of the charging curve, but it can fail for large 
values of the Debye screening length, since in this case 
further-nearest neighbor interactions become important. 
For large long-range interactions can be important 
although they are much weaker than the interactions be- 
tween neighboring monomers, as shown in the following 
numerical examples. 



1. Results for radially symmetric charge distribution 

In Fig. 4(a) the NN prediction (dashed line, exact 
solution) is compared with an enumeration using the 
full long range interaction (circles). The Debye length 
is = 100 A and interactions between monomers 

are calculated using Eq. ^ with b = 0, d = 2.5 A, 
£w/£d = 80/3 and a monomer separation a — 3.5 A . 
In all calculations with dielectric discontinuity we use a 
monomer separation a = 3.5 A rather than 2.5 A, which 
corresponds to a somewhat smaller fraction of dissociable 
groups. 

The exact solution of the NN model in Fig. 4(a) (bro- 
ken line) deviates significantly from the enumeration re- 
sults with the full range of interactions included (cir- 
cles). In contrast, our non- uniform mean- field approach 
with symmetry breaking (solid line in Fig. 4(a)) is semi- 
quantitatively correct. The success of our generalized 
mean- field approximation is one of the main results in the 
first part of this work. Note that the main difference with 
respect to enumeration is that the mean-field approxima- 
tion over-estimates the effects of symmetry breaking, as 
seen from the exaggerated size of the plateau region. 

We also present in Fig. 4(a) a comparison between 
the NN exact solution (dashed line) and an enumera- 
tion taking only nearest-neighbor interactions into ac- 
count (crosses). This is done in order to test finite size 
effects in the exact enumeration. The enumeration and 
exact solution yield almost identical results, demonstrat- 
ing that an enumeration with = 20 is typically very 
accurate for a single polymer. We note that periodic 



boundary conditions are essential in order to obtain this 
level of accuracy in enumeration, when long range in- 
teractions are included (compare Fig. lb). As another 
test for the enumeration procedure we increased N to 30 
for several different choices of vdb_{z). In all these cases 
deviations from results with A'^ = 20 were insignificant. 

Figure 4(b) shows a comparison between the NN exact 
solution (broken line) and a mean-field calculation with 
two sublattices (solid line), taking only the NN interac- 
tion into account, for n^^ = 100 A . The NN-mcan-field 
approximation shows an effect similar to the curve in 
Fig 4(a): two artificial discontinuities in the derivative of 
(s) are seen, corresponding to two erroneous second-order 
transitions. In addition, the plateau is more pronounced 
compared to the exact solution. Nevertheless, the overall 
prediction for (s) as function of pH is quite accurate. 

Finally, the dependence on the Debye screening length 
is investigated in Fig. 4(c), using our non-uniform 
mean-field approach. The parameters J and K of the 
sublatticc interactions are calculated for each value of k 
using Eqs. H12|) and H23|) . For all the four values of 
that are shown, = 3, 10, 30 and 100 A a pronounced 
plateau is visible. The cusps that are present in all four 
cases are artifacts due to the mean-field approach. 



2. Dependence on the position of charged groups 

In most polyacids and polybases the charged units are 
located in side groups, rather than being close to the 
polymer axis. This raises the question whether an in- 
teraction much larger than the usual Debye-Hiickel in- 
teraction will occur even if the charges are displaced 
from the axis. Within the simple cylindrical model pre- 
sented above, this question can be addressed by calculat- 
ing the electrostatic interaction between two monomers 
as a function of b. Such a calculation is shown if Fig. 5, 
for four different values of the screening length k^^. The 
monomers are separated by a distance of a = 3.5 A, while 
the polymer radius is taken as d = 2.5 A, as in Fig 4. In 
all four cases a very large decrease of wdh occurs with in- 
crease of b toward the cylinder boundary, b — r. For large 
r this result is not surprising, since the cylinder becomes 
similar to a planar interface, separating an aqueous ionic 
solution and a low dielectric medium. Near such an in- 
terface the electrostatic interaction is equal to twice the 
screened electrostatic interaction in water . For large 
cylinder radius d we have checked that Eq. 

indeed yields this result. In Fig. 5 the cylinder radius is 
not large compared to the Debye length and the interac- 
tion close to the cylinder boundary is even smaller than 
in the planar limit. The screened interaction in water is 
shown for comparison using dashed lines. 

The above analysis demonstrates that actual electro- 
static interactions between near-by monomers depend 
strongly on the spatial organization of the PE. These 
interactions probably cannot be estimated reliably using 
simplified models such as the cylindrical one presented 
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FIG. 4: Average degree of dissociation of a polymer with strong interactions between close by monomers, characterized by a 
plateau near (s) = 1/2. The interactions between monomers are calculated using the cylindrical model shown in Fig. 3, with 
£d — 3, Ew = 80, d = 2.5 A, a — 3.5 A, fe = and a Debye length — 100 A. (a) Mean-field results with two sublattices 
(solid line) are compared with an enumeration over all configurations with N — 20 and periodic boundary conditions (circle 
symbols). The mean-field results show two cusps, where transitions occur between a solution with no symmetry breaking and 
one with symmetry breaking (the latter occurs for intermediate pH). Results are also compared with an enumeration taking 
only nearest neighbor (NN) interactions into account (crosses) and the exact solution with NN interactions, calculated using 
the transfer matrix method (dashed line), (b) Mean-field calculation with two sublattices taking only NN neighbor interactions 
into account (solid line), compared with the exact solution with NN interactions (dashed line), (c) Dissociation as function of 
PH for four different values of the Debye screening length: 3, 10, 30 and 100 A, calculated using the mean-field approximation 
with two sublattices. All parameters other than k are as in part (a). The interaction parameters J,K used in the mean-field 
approximation are equal to 0.21,4.36 (k'^=3A); 0.46,4.86 (^"^ = 10 A); 0.87,5.30 (k"^=30A); and 1.37,5.84 («;-^ = 100A). 



above. The detailed polymer structure, as well as other 
effects such as the discreteness of the solvent, must be 
taken into account. 



D. Further discussion of the two-sublattice 
approximation 

Although the comparison with exact enumeration 
demonstrates that a two-sublattice model is useful, one 
may ask to what extent the separation into two sublat- 



tices has a physical significance. In order to discuss this 
question we note that within the plateau region of the 
titration curve there is typically a strong anti-correlation 
between even and odd monomers. In order to under- 
stand this anti-correlation one may think of the ground 
state of the Hamiltonian in Eq. ((BJ, concentrating first 
on the special case // = 0. When the interaction be- 
tween monomers favors opposite dissociation values, the 
ground state is typically a periodic array of alternating 
values in the even and odd positions, s = +1 and s = —1. 
The long range order of the ground state is not preserved 




Ql 1 1 1 1 1 

0.5 1 1.5 2 2.5 

b[A] 

FIG. 5: Dependence of the electrostatic interaction between 
two monomers on their distance 6 from the polymer axis, 
within the cylindrical model of Fig. 3. The distance between 
the monomers is a = 3.5 A and the other model parameters 
are d = 2.5 A, £d = 3 and e„ = 80. Results are shown for four 
values of the Debye length, =3, 10, 30, and 100 A. The 
electrostatic interaction in an aqueous ionic solution is shown 
for comparison using dashed lines. 



within the exact theory at any finite temperature, due to 
the entropy associated with domain boundaries in a one 
dimensional system ^TTj. Nevertheless, at a certain range 
of pH values around /i = we may expect a staggered 
correlation function with strong anti-correlation between 
even and odd monomers. 

As an example, Fig. 6 shows the correlation function, 
(sjSo) — (si) (so); for a PE having the same parameters 
as in Fig. 4(a), calculated by exact enumeration over all 
states of a PE of length N = 30 and using periodic 
boundary conditions. Results are shown for three dif- 
ferent pH values: in the top plot pH = 7, corresponding 
to /i = 0. The correlation function has a staggered form 
which persists over the full length of the PE. Nevertheless 
it is clear that there is no true long range order because 
the (anti) correlation decreases slightly with monomer 
separation. In the middle plot, where the pH is equal to 
9, the correlation has a shorter range, persisting only up 
to a distance of about 8 sites from the center monomer. 
Note that a pH value of 9 is approximately at the right 
edge of the plateau region seen in Fig. 4(a). With fur- 
ther increase of pH the correlation length continues to 
decrease and at pH = 13.5 (bottom plot) there is almost 
no correlation even between adjacent monomers. This pH 
value is close to the transition point that is found in the 
two-sublattice model, beyond which there is no symme- 
try breaking between the two sublattices. In summary. 
Fig. 6 demonstrates that the two-sublattice approxima- 
tion captures an essential physical property of the disso- 
ciation pattern that is absent in the uniform mean field 
approach, namely a strong anti-correlation between even 
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FIG. 6: Correlation function between the dissociation of a 
monomer and that of its neighbors, (siSo) — (si) (so), calcu- 
lated from an exact enumeration over all configurations of a 
PE having A'^ = 30 monomers, with periodic boundary condi- 
tions. All physical parameters are as in Fig. 4(a). The pH is 
equal to 7 in (a), corresponding to /i = 0, to 9 in (b), and to 
13.5 in (c). 



and odd sites. 



In principle, a periodicity other than two may be in- 
cluded in the formulation of the mean field equations, and 
could lead, for certain parameters, to a lower free energy 
than the two-fold periodicity. In such cases a plateau 
would be expected in the titration curve at an average de- 
gree of dissociation other than one half. Comparison with 
the exact enumeration and with typical experimental re- 
sults indicates that such additional symmetry breaking 
into structures with more than two sublattices does not 
occur within the physical parameters considered in this 
work. 



In the second part of this work, where we will look 
at the interaction between two weak polyelectrolytes, we 
will employ simple Debye-Hiickel interactions, as well as 
the interaction within a cylindrical dielectric cavity with 
6 = 0, which constitute the two extreme cases. In both 
cases we expect a mean-field approach with two sublat- 
tices to be adequate in order to predict the average charg- 
ing, as was demonstrated in the preceding discussion. 
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FIG. 7: Schematic illustration of a model describing interac- 
tion of a polyacid with a parallel polybase, separated by a 
distance D. The distance between charged groups in both of 
the PEs is a. For simplicity the charged groups are facing 
each other. 



III. INTERACTION BETWEEN POLYACID 
AND POLYBASE 

A. Uniform mean-field approach 

The model we consider is shown schematically in Fig. 7. 
A polyacid (left) and polybase (right) are aligned paral- 
lel to each other and separated by a distance D. For 
simplicity we assume that the distance between charged 
groups (denoted by a) is identical in the two polymers 
and that the charge lattice are in phase with each other 
in the two polymers, as shown in the figure. We would 
like to calculate the average charge on the two polymers 
and the free energy as a function of D. 

We consider first the case of a uniform dielectric con- 
stant, Eq. Q), and also assume that for each polymer 
a uniform mean-field theory (with no symmetry break- 
ing) is adequate. As was shown in the calculations for 
a single polymer, the latter assumption is justified for 
monovalent monomers having a nearest neighbor separa- 
tion a > 2.5 A. 

It is convenient to define for the polyacid 



1 



and for the polybase 



Sb 



1- h 



(26) 



(27) 



where Sa and are zero for an uncharged monomer and 
one for a charged (dissociated) one. With these defi- 
nitions both Sa and Sb increase with pH. The mean-field 
equations are found in a similar way as in the single poly- 
mer case, and are given by 



ha = fj-a + J (Sa)o + K 

hb = fj-b + J {sb)o + K (Sa)^ 



{sa)o = -ta.nh{ha) ; {sb)o = -tanh(/ih). 
The coefficients in these equations are given by: 



(29) 



J = 



2.303 

2~ 
2.303 



(pH - pKa) + Afl 
(pH - pKb) - A/i 



(30) 



where 



(31) 

Equations (|2HJl-l|2ni) are very similar to Eqs. ((Tm) - ((TT|l . 

with a number of important differences. First, the sub- 
scripts a and b do not represent two sublattices but 
instead distinguish between the polyacid and polybase. 
Another difference is that the chemical potentials jla and 
jib are usually not equal to each other. Most importantly, 
J is almost always larger than if, whereas for the two 
sublattice case J is smaller than K. It is easy to show 
that for J > if Eqs. have a single solution. 



1 . Results 

The electrostatic interaction between the polyacid and 
polybase increases dissociation in both polymers (in con- 
trast to the interactions within each PE, which inhibits 
charged groups from dissociating). Figure 8 shows the 
degree of charging of a polyacid and polybase having 
a — 2.5 A, as function of pH and for three different values 
of D. When the polymers are sufficiently far away from 
each other their dissociation curves are identical to those 
of a single polymer. For smaller separation the average 
charge increases. An important case occurs when the pH 
is tuned such that 



(pH - pKa) = - (pH - pKb) 



(32) 



For example, with the parameters used in Fig. 8, pKa = 
4 (similar to poly-acrylic-acid) and pKb = 10 (similar to 
poly-vinyl-amin), this equality holds at pH — 7. In this 
case, one has /it, — —jla, as seen from Eq. H30I) . and the 
solution of Eqs. H28|l - H29|) has the properties ha = —hb 
and (sq)q = — (sb)o- Using the definitions in Eqs. if^ - 
(|27|l . the average charging degrees of the polyacid and 
polybase are then equal to each other, and the value of 
ha is found from the single transcendental equation: 



(28) 



ha = fia + {K - J) tanh(/ia) 



(33) 
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FIG. 8: Average dissociation degree of a polyacid (increasing 
lines) interacting with a polybase (decreasing lines), as func- 
tion of pH. The distance between charges is a = 2.5 A; pKa 
= 4, pKb = 10 and the Debye length is = 30 A. Results 
are shown for three values of the inter-polymer separation: 
7? = 100 A (solid lines), 10 A (dashed lines) and 5 A (dash- 
dot lines) . The symmetric case where pH — pKa = pKb — pH 
occurs at pH = 7. 



In the following examples we restrict ourselves to the 
symmetric case described by Eqs. (|32f) and Figure 
9(a) shows the average degree of dissociation as func- 
tion of the polymer separation D (identical for the poly- 
acid and polybase). Results are shown for a — 2.5 A, 
pH — pKa = pKb — pH = 3 and for four different val- 
ues of the Debye length, ranging between 3 A and 100 A. 
When D is large compared to the polymers do not 
interact, and their average charge is equal to its value 
in an isolated polymer (compare with Fig. 2 at pH-pKa 
= 3). This value depends strongly on k~^. At separa- 
tions D of order k^^ and smaller, the average charging 
increases with decrease of D and approaches unity (full 
dissociation) at contact. 

We turn to the free energy of the two interacting PEs, 
shown in Fig. 9(b). In this figure the free energy F is 
divided by N, the number of monomers in each PE. A 
distinctive feature in this figure is that F is almost inde- 
pendent on at small separations. In contrast to this 
short separation behavior, F depends strongly on k^^ at 
large separations. In order to understand these two be- 
haviors we consider each one of the two limits separately: 



2. Small PE distances, D <^ n ^ 

At short separations the average degree of dissociation 
saturates and is independent on as can be seen in 
Fig. 9(a). The electrostatic interaction energy of the two 
polymers also becomes nearly independent on k^^, as 
can be understood from the following argument. Con- 
sider the two polymers as uniformly, oppositely charged 



and parallel lines. In the limit of no screening, k = 0, 
the electrostatic energy is dominated by interactions at 
distances of order D and smaller. At distances larger 
than D opposing positive and negative charges in the two 
polymers can be regarded as dipoles and the electrostatic 
interaction between them decays as 1/z^ where z is their 
distance, measured parallel to the polymers. As long as 

^ D screening affects only these dipole-dipole inter- 
actions, but not the main electrostatic contribution com- 
ing from interactions at distances smaller than D. The 
independence of both (s) and the electrostatic energy on 

leads to the behavior seen at these small separations. 



3. Large PE distances, D ^ k ^ 

At large separations the free energy approaches the 
sum of free energies of the two isolated polymers. More 
precisely, the average degree of dissociation on the two 
polymers approaches a constant and the free energy can 
be approximated as follows: 



F 

N 



2Fo 

N 



21, 



(^^yKoi^D) (34) 



where Fq and sq are the free energy and average degree of 
dissociation of a single, isolated PE, respectively. These 
constants are unrelated to the interaction between the 
two PEs but depend strongly on k. The modified Bessel 
function Ko{kD) characterizes the electrostatic interac- 
tion between two parallel and uniformly charged rods: 



Ko{kD) 



dz 



exp (— K\/ + D^) 

+ 



' exp {-Vu^ + K^D'^) 
au = — - (35) 

x/u2 + «;2£,2 



Deviations from the asymptotic form (|34|l are expected 
to occur only when the average degree of dissociation 
deviates from sq. This happens approximately when 
D < as can be seen in Fig. 10, where the differ- 
ence between F/N and Eq. (|34l) is plotted for the same 
four values of as in Fig. 9. 



B. Non-uniform mean-field approach 

In the case of stronger interactions between monomers 
within each polymer, the dissociation curve of a single 
PE is characterized by a plateau. In this case we expect 
a symmetry breaking transition with two sublattices on 
each PE. In order to deal with this case the mean- field 
Hamiltonian can be generalized to account for sublattices 
on each one of the two PEs: 



■Hi 



(36) 
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FIG. 9: Average degree of dissociation (a) and free energy per monomer, F/N (b) as function of the distance D between a 
polyacid and polybase, with pH — pKa = pKb — pH — 3. Simple Debye-Hiickel interactions are used with = 3 A (dashed 
line), 10 A (solid line), 30 A (dotted line) and 100 A (dash-dot line). All other parameters are as in Fig. 8. The arrows on the 
right hand side of (a) show the value of (s) for an isolated PE. 
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FIG. 10: The difference AF between the free energy F of a 
polyacid interacting with a polybase, and the asymptotic form 
of Eq. 1341 . Four values of are shown. These values, and 
all other parameters and notations are as in Fig. 9. For each 
value of the values of Fo and so in Eq. 13411 are equal 
to the free energy and average degree of dissociation of an 
isolated PE, respectively. 



Eq. l|?7|l are equal to: 



i#0 

^ = ^E^DH[2(* + l)a] 

i 
i 

h = ^^^^DH [V(2i + l)2a2 + £)2 



(38) 



and p}' are given by Eqs. (|30|l and (|31|l . The four 
equations for hf typically have multiple solutions; For ex- 
ample, if symmetry breaking occurs on both polymers the 
number of solutions is 9. In the limit of non-interacting 
polymers, Iq = Ii = 0, four of these solutions are equiv- 
alent minima of the free energy, related to each other by 
exchange of the two sublattices on one or both of the 
polymers. Interactions between the polymers break the 
symmetry of exchanging only the sublattices in one of the 
polymers, and there are two (equivalent) global minima 
of the free energy. 



The mean-field equations and free energy are found using 
the Gibbs variational principle, in a similar fashion as for 
the single polymer case. For example, the equation for 
/ig is: 

/ig = A'' -I- J (So^o + K (h)^, + Jo (So)o + h (Si)l (37) 

where {si)^ = — tanh/if . Similar equations are obtained 
for hi, Hq and h\. The coefficients J,K,Io and Ii in 



1. Results 

As a concrete example we consider again the model 
shown in Fig. 3, which accounts for a low dielectric con- 
stant of the polymer backbone. Parameters are similar 
to Fig. 4, & = 0, a = 3.5A, d = 2.5 A and k'^ = 100 A. 
The pH, pKa and pKb values are chosen such that 
pH — pKa = pKb — pH = 3. The coefficients J and 
K are set as in Fig. 4 with = 100 A, J = 1.4 

and K — 5.8. For the coefficients Iq and Ii we use 
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FIG. 11: Average degree of dissociation of interacting poly- 
acid and polybase, as function of their distance. The inter- 
actions between monomers within each PE are calculated as- 
suming a low dielectric backbone, as in Fig. 3. The param- 
eters of the model are as in Fig. 4: a = 3.5 A, d = 2.5 A, 
6 = and = 100 A, and pH - pKa = pKb - pH = 3. The 
inset shows the free energy per monomer, F/N, as function 
of D (solid line). The dashed line shows the approximation 
of Eq. 1)34^ with sq = 1/2 and Fq set to match the value of F 
at large D. 



Eqs. 1)38(1 with the screened Debye-Hiickel interaction in 
water, Eq. Q). Both of these choices are approxima- 
tions which become inaccurate when the polymers are 
very close to each other, since the electrostatic Green's 
function should then be evaluated in the presence of two 
dielectric cylinders. However we expect our results to be 
quahtatively correct as will be further discussed below. 

Figure 1 1 shows the average charging of the two poly- 
mers as function of their separation D. Due to the 
plateau in the dissociation curve of each polymer, (s) is 
close to 1/2 in most of the separation range. A sharp 
increase in (s) is found at close separations of a few 
Angstroms. Note that this range, where interactions be- 
tween the polymers affect the average charging, is much 
smaller than the Debye length, = 100 A. 

The inset shows the free energy (solid line) as function 
of D. For comparison we show (by a broken line) the 
approximation of Eq. H34|l . with sq equal to 1/2 and Fq 
matching the value of F at large D. We note that the 
constant Fq used in Fig. 11 is not exactly the free energy 
of an isolated PE. At very large separations, when the 
PEs are truly isolated from each other, the average degree 
of dissociation is somewhat smaller than 0.5 (see Fig. 4 
at pH-pKa = 3). At these separations the asymptotic 
form of Eq. H34|l applies with Fq and sq set according to 
the properties of the isolated PEs, while the form shown 
in Fig. 11 matches the free energy at intermediate values 
of D where (s) « 0.5. The free energy F deviates from 
Eq. (|34|) only at very large separations, and at very small 
separations where (s) is larger than 1/2, as seen in the 



inset of Fig. 11. 

In summary, the symmetry-broken solution is stable 
for a wide range of distances and gives way to a sym- 
metric solution only at very small distances. With the 
dielectric discontinuity on both polymers taken properly 
into account (affecting Iq and Ii as well as K and J) we 
can expect a stronger interaction between the polymers 
for small D. This will lead to an increase of the average 
charging at somewhat larger values of D than in Fig. 11. 



IV. SUMMARY 

The main result in the first part of this work concerns a 
generalization of the standard mean-field theory of charge 
regulation in weak PEs. The polymer is divided into 
two sublattices, allowing explicitly for correlations be- 
tween these sublattices to be taken into account. Similar 
models have been studied in the past in the context of 
Ising-like models. In the Ising model, interactions are 
usually assumed to be short-ranged. If only interactions 
between neighboring monomers are considered, the par- 
tition function can be calculated exactly. For PEs the 
main advantage of using a mean-field approximation is 
that it allows long-range electrostatic interactions to be 
taken into account. Simultaneously, one expects mean- 
field methods to gain in accuracy as the range of inter- 
actions increases. Our main result is that a mean-field 
approach with separation into two sublattices is adequate 
within a wide range of model parameters. In particular it 
succeeds in the case of large inter-monomer interactions, 
where a uniform mean- field theory fails, while also taking 
into account long range interactions, which may still play 
an important role. 

A motivation for the use of a nearest-neighbor approx- 
imation was recently suggested in Ref. *?! . It was pointed 
out that a low dielectric constant of the polymer back- 
bone can lead to strong enhancement of the coupling be- 
tween close-by monomers. We show that even within 
the model of Ref. Q, the nearest neighbor approxima- 
tion needs improvement for large values of the Debye 
length, because of the contribution of interactions be- 
tween non-neighboring monomers. On the other hand, 
a mean-field approximation with two sublattices is semi- 
quantitatively accurate. We also demonstrate that effects 
due to the dielectric discontinuity between the PE inte- 
rior and the aqueous solvent depend sensitively on the 
location of the charged groups within the low-dielectric 
cavity; this is quite relevant, since for most experimen- 
tal PE architectures, the charged groups are not located 
centrally but are displaced towards the aqueous interface. 

The linearized Dcbyc-Hiickel theory is used in this 
work to evaluate the interaction between monomers. 
This, of course, is only an approximation, whereas in 
principle the full non-linear response of the ionic solu- 
tion must be taken into account. Use of Debye-Hiickel 
interactions is justified as long as the electrostatic po- 
tential is small compared to the thermal energy. Hence 
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this approximation is probably reasonably accurate in the 
plateau region, where the average charge along the poly- 
mer is small. Far away from the plateau, and for highly 
charged PEs, one needs to go beyond Debye-Hiickel the- 
ory, using the non-linear Poisson-Boltzmann theory. The 
great advantage of using pairwise interactions is that they 
allow tractable, analytic solutions to be obtained. In con- 
trast, the nonlinear distribution of ions cannot be calcu- 
lated analytically even near a uniformly charged cylinder 
immersed in a salt solution, let alone an inhomogeneously 
charged PE. In light of this situation we believe that the 
results presented in this work provide a useful qualitative 
treatment of charge regulation even for the case of highly 
charged PEs, although they may be quantitatively mod- 
ified by charge renormalization due to nonlinear effects 
close to the PE [l^. 

In the second part of this work we studied the in- 
teraction between a polyacid and a polybase, using a 
mean-field approximation. This interaction leads to an 
increase of the average charging in both polymers as 
they approach each other. In addition, the interaction 
energy between a weak polyacid and a weak polybase 
is stronger than expected in the absence of distance- 
dependent charge regulation. This means that the ef- 
fect of charge regulation may be important for the build- 
ing of stable multilayers, since it allows to simultane- 
ously decrease the repulsion between similarly charged 
weak PEs (since the charge regulation in this case de- 
creases the charge strength) but also leads to strongly 



bound polyacid-polybase pairs. For close-by polymers 
the electrostatic energy is dominated by interactions be- 
tween neighboring monomers, and the free energy de- 
pends only weakly on the Debye screening length k"^. 
On the other hand at large separation between the poly- 
mers both the average degree of dissociation and the free 
energy typically vary strongly with k^^ . The characteris- 
tic distance where interactions between the polymers can 
affect their degree of dissociation is the Debye screening 
length. However, when there is a plateau in the charge vs. 
pH curve of a single PE, the degree of dissociation may 
remain close to 1/2 even at small separations compared 
to . In these cases a sharp increase in the average 
dissociation can occur close to contact. 

The increase of charging when weak polyacids and 
weak polybases come into contact could in principle 
be observed using infra-red spectroscopy in multilayers. 
However, one has to keep in mind that in such highly 
concentrated systems the oppositely charged groups will 
get very close to each other and form salt bridges. 
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